Introduction

The absence of translational periodicity and symmetry, and the rich structural complexity make it difficult to understand the order within disorder1,2 in disordered materials. The advent of advanced instrumentation and measurement protocols makes it feasible to use quantum beam diffraction (X-ray diffraction (XRD) and neutron diffraction (ND)) techniques to reveal the structure of disordered materials at synchrotron and/or neutron sources3,4,5,6.

Amorphous (a)-silicon (Si) and silica (SiO2) are the most typical and important disordered materials in both fundamental and technological research studies. In particular, SiO2 is a canonical amorphous-former, whereas it is possible to synthesize a-Si only in a thin film owing to its poor amorphous-forming ability. The short-range structural unit of these amorphous materials is a tetrahedron, SiSi4 in a-Si and SiO4 in a-SiO2, and the formation of their networks is governed by the corner-sharing of tetrahedra. This corner-sharing motif is within Zachariasen’s classification7 for the glass formation of oxide materials, but elemental amorphous materials do not follow the rule, because it is impossible to obtain bulk amorphous silicon as mentioned above.

The structures of these materials have been studied by diffraction techniques. Laaziri et al. reported high-quality X-ray diffraction data for a-Si with a high-real space resolution for precisely determining the coordination number8. The structure of a-SiO2 has been widely studied by both X-ray9,10 and neutron diffraction4,11,12. The most important feature in the diffraction data of a-SiO2 is that the first sharp diffraction peak (FSDP)11,12,13,14 is observed in both X-ray and neutron diffraction data, whereas the second diffraction peak, the so-called principal peak (PP)12,14, can be observed in only the neutron diffraction data because this peak reflects the packing of oxygen atoms12,15, which is sensitive to neutrons. The origin of the FSDP has been discussed for a long time. The FSDP was first discussed in 197616, although it seems that the name “FSDP” was first used by Phillips in 198117. The interpretation of diffraction peaks including the FSDP was attempted in the 1980s13,17, as discussed in details in several papers18,19,20,21,22. It is known that the FSDP of a-SiO2 is related to the formation of the random network model proposed by Zachariasen7, which was extended to silicate glasses by Greaves23 and recently revised by Mei et al., as illustrated in Ref. 21. It was demonstrated that intermediate-range ordering arises from the periodicity of boundaries between successive cages in the network formed by the connection of regular SiO4 tetrahedra with shared oxygen atoms at the corners associated with the formation of a ring structure and a large cavity21,22. The second maximum, PP, reflects the size of the local-network-forming motif, whereas the FSDP indicates the arrangement of these motifs in an intermediate range according to Zeidler and Salmon24. Another interpretation of the FSDP has recently been proposed by Shi and Tanaka, who discussed local tetrahedral ordering in covalent liquids and glasses25 and they concluded that a-Si has an FSDP. Moreover, an FSDP was found in metallic glasses26, although Price et al. implied that some diffraction peaks from amorphous alloys are not FSDPs13.

The investigation of the behavior of the FSDP and PP under high pressures is important to understand the nature of intermediate-range ordering in disordered materials. Figure 1a shows in situ neutron structure factors, S(Q), for a-SiO2 under high pressures reported by Zeidler et al.27. It was found that the FSDP diminishes with increasing pressure, simultaneously with a peak shift to a higher Q. On the other hand, the PP becomes sharp with increasing pressure, suggesting that the oxygen packing fraction increases12,15 with the decrease of cavity volume under high pressures. This is an important benchmark for the modification of the intermediate-range ordering of a-SiO2 under cold compression. Onodera et al. recently reported the unusual behaviour of the FSDP after hot compression2; in the X-ray S(Q) shown in Fig. 1b, they have observed the evolution of FSDP at 7.7 GPa at a temperature higher than 400 °C. These diffraction data suggest that FSDP is very sensitive to pressure and temperature, while the PP position is insensitive to the density change.

Figure 1
figure 1

(a) In situ neutron structure factors, S(Q), for a-SiO2 under cold compression reported by Zeidler et al.27 (b) X-ray structure factors, S(Q), for a-SiO2 after hot compression reported by Onodera et al.2

In this article, we apply several topological techniques to analyze the ring size and cavity distributions, and tetrahedral order of crystalline, amorphous, and liquid Si and SiO2 to reveal the network topology2,28,29,30 for understanding the diffraction peaks in disordered materials with a special focus on the FSDP and PP.

Methods

Structure modeling

Atomistic models of liquid (l)-SiO2 and a-SiO2 were obtained by classical molecular dynamics (MD) simulation and MD–reverse Monte Carlo (RMC) modeling29, respectively. The atomistic models for l- and a-SiO2 were generated by combined classical MD simulation–RMC modeling. The MD model for l-SiO2 was adopted from the literature31. The MD simulation for a-SiO2 was performed using the large-scale atomic/molecular massively parallel simulator (LAMMPS) package32 within the NVT ensemble. The interactions were described by pair potentials with short-range Born–Mayer repulsive and long-range Coulomb terms, i.e.,

$${\phi }_{ij}\left(r\right)={B}_{ij}{\text{exp}}\left(-\frac{r}{{\rho }_{ij}}\right)+\frac{{e}^{2}}{4\pi {\varepsilon }_{0}}\frac{{Z}_{i}{Z}_{j}}{r},$$
(1)

where r is the distance between atoms i and j, Bij and ρij define the magnitude and softness of the Born–Mayer terms, respectively, Zi is the effective charge on atom i (ZSi = 2.4, ZO = -1.2), e is the elementary charge, and ε0 is the permittivity of vacuum. The Bij values were 21.39 × 10−16 J (Si–O), 0.6246 × 10−16 J (O–O) or zero (Si–Si); the ρij values were 0.174 Å (Si–O), 0.362 Å (O–O) or zero (Si–Si). As initial configuration, 3000 (Si, 1000; O, 2000) atoms were randomly distributed in a cubic cell with a side length of 35.66 Å. The cell had a number density of 0.06615 Å-3. The simulation used periodic boundary conditions, and the long-range Coulomb interactions were treated by using the Ewald summation. A time step of 1 fs was used in the Verlet algorithm. The simulation temperature was maintained at 4000 K for 20,000 time steps, then the temperature was reduced to 300 K over 200,000 time steps. Finally, the system was equilibrated at 300 K for 50,000 time steps. A Nosé-Hoover thermostat was employed to control the temperature. After MD simulations, the configurations obtained for l- and a-SiO2 were refined by RMC modeling with constraints on the coordination numbers and on the O–Si–O bond angle distribution, in order to prevent the formation of an unfavorable disordered structure.

Atomistic models of l-Si and a-Si were obtained by RMC modeling29 and combined classical MD simulation–RMC modeling, respectively.

The RMC model for l-Si (1770 K, 5000 particles) was obtained by the RMC + + code33 and based on the X-ray S(Q) for l-Si29. The number density was 0.055 Å−3, that is consistent with a bulk mass density of 2.57 g cm−334.

The model for a-Si was obtained by MD simulation, followed by RMC refinement. The MD simulation was performed using the LAMMPS code. The modified Tersoff potential 35 based on the three-body Tersoff potential 36,37 was used for describing interatomic interactions. In the modified Tersoff potential, the total energy E is written as follows:

$$E=\frac{1}{2}\sum_{i}\sum_{j\ne i}{V}_{ij}({r}_{ij}),$$
(2)
$${V}_{ij}\left({r}_{ij}\right)={f}_{C}\left({r}_{ij}\right)\left[A{\mathrm{exp}}\left(-{\lambda }_{1}{r}_{ij}\right)-{b}_{ij}B{\mathrm{exp}} \quad \left(-{\lambda }_{2}{r}_{ij}\right)\right],$$
(3)
$${f}_{C}\left(r\right)=\left\{\begin{array}{c}1, if \quad r \le {R}_{1}\\ \frac{1}{2}\left[1+\mathrm{cos}\left(\pi \frac{r-{R}_{1}}{{R}_{2}-{R}_{1}}\right)\right]\\ 0, if \quad r \ge {R}_{2}\end{array}\right., if \quad {R}_{1}<r<{R}_{2},$$
(4)
$${b}_{ij}={\left(1+{\xi }_{ij}^{n}\right)}^{-\frac{1}{2n}},$$
(5)
$${\xi }_{ij}=\sum_{k\ne i,j}{f}_{C}({r}_{ij}){g}_{\text{mod}}\left[{\theta }_{ijk}\right]{\text{exp}}\left[\alpha {\left({r}_{ij}-{r}_{ik}\right)}^{m}\right],$$
(6)
$${g}_{\text{mod}}\left(\theta \right)={c}_{1}+{g}_{o}\left(\theta \right){g}_{a}\left(\theta \right),$$
(7)
$${g}_{o}\left(\theta \right)=\frac{{c}_{2}{\left(h-\mathrm{cos} \theta \right)}^{2}}{{c}_{3}+{\left(h-\mathrm{cos} \theta \right)}^{2}},$$
(8)
$${g}_{a}\left(\theta \right)={1+c}_{4}{\text{exp}}\left[-{c}_{5}{\left(h-\mathrm{cos} \theta \right)}^{2}\right],$$
(9)

where E is decomposed into bond energies Vij between atoms i and j. rij is the distance between atoms i and j, θijk is the angle confined by the bonds between ij and ik. The function Aexp(-λ1rij) and bijBexp(− λ2rij) represent a repulsive and an attractive term, respectively. The extra term fC is merely a smooth cutoff function, to limit the range of the potential. bij represents a measure of the bond order. The term 1 + ξij corresponds to the coordination number of atom i, where ij bond and other bonds are counted as 1 and ξij respectively. The term gmod(θ) is the modified angular-dependent term. The potential parameters, A, B, λ1, λ2, n, α, m, c1, c2, c3, c4, c5, h, R1, and R2, are given in Table 1. The simulation box was a cubic cell with a side length of 50.00 Å. In the MD simulation, 6256 Si atoms were placed in the box as an initial configuration with in the NVT ensemble. A time step of 1 fs was used in the Verlet algorithm. The atomic configuration was initialized at random and the system was equilibrated at 3000 K for 500,000 steps. Then it was cooled to 300 K during 5,000,000 steps and annealed at 300 K for 500,000 steps. After the MD simulation, the configuration obtained was refined by RMC simulation with constraints on the coordination number and the Si–Si-Si bond angle distribution, in order to avoid formation of unfavorable disordered structures. The model reproduces the experimental X-ray S(Q).

Table 1 Potential parameters of the modified Tersoff potential.

For l-Si, additional principles (FP)MD simulations within the framework of density functional theory were performed for a 64-atom Si supercell with periodic boundary conditions. The calculations were performed using the projector-augmented wave method38 and the generalized-gradient approximation with the exchange correlation functional of Perdew, Burke, and Ernzerhof39. The electronic wave functions were expanded in a plane-wave basis with an energy cutoff of 400 eV at the Γ point in the Brillouin zone. Atomic configurations selected from the FPMD trajectory generated using the Vienna Ab initio Simulation Package40 were used to validate the atomistic model of l-Si.

Ring size distribution analysis

The ring size distribution was calculated using the R.I.N.G.S. code41,42 on the basis of the King43 and primitive44,45 criteria. The first coordination distances were set to 2.5 Å and 3.0 Å for a-Si and l-Si, and to 1.9 Å and 2.4 Å for a-SiO2 and l-SiO2, respectively.

Cavity volume calculation

Cavity volume analyses were performed using the pyMolDyn code46. The code can calculate three different types of cavity, domain, center-based (Voronoi), and surface-based cavities. We calculated surface cavity volumes with a cut off distance rc of 2.5 Å.

Tetrahedral order parameter q

The tetrahedral order parameter for the Si-centered hyper-tetrahedra and SiO4 tetrahedra are defined as47

$$q\equiv 1-\frac{3}{8}\sum_{i=1}^{3}\sum_{k=i+1}^{4}{\left(\frac{1}{3}+{\mathrm{cos}} {\theta }_{ijk}\right)}^{2},$$
(10)

where θijk is the angle formed between the central Si atom j and its neighboring Si or O atoms i and k. This parameter was designed to give a value of unity for a regular tetrahedron and a mean value of zero for a perfect gas.

Results and discussion

The X-ray structure factor, S(Q), for a-Si8 is shown in Fig. 2a together with that of l-Si (1770 K)29. We observe Q2 (PP) and Q3 at QrA-X ~ 5 and 8.5, respectively, and no Q1(FSDP) is observed. Note that scattering vector Q is scaled by multiplying by rA-X (first coordination distance)6,11,12,13,14,15,18,22,24,29,48 obtained by a Fourier transform of S(Q). S(Q) for l-Si is different from that for a-Si owing to the increased coordination number, NSi-Si, from approximately 3.9 (amorphous)8 to 5.7 (liquid)29 associated with the significant density increase from 2.30 g cm−3 (amorphous) to 2.57 g cm−3 (liquid). This behavior was also confirmed in previous FPMD studies49,50 and is consistent with the fact that a-Si is a semiconductor and l-Si is a metallic liquid. X-ray29,30 and neutron30 structure factors, S(Q), for a-SiO2 are shown in Fig. 2b. As mentioned in the previous section, Q1 (FSDP) is observed for a-SiO2 in both the X-ray and neutron S(Q), but Q2 (PP) is visible only in the neutron S(Q), because it reflects the packing fraction of oxygen atoms12,15 since neutrons are sensitive to O–O correlations but X-rays are more sensitive to Si–Si correlations. It is worth mentioning that the X-ray S(Q) of l-SiO2 (2019 K) is identical to that of g-SiO2. Q1 (FSDP) in the X-ray S(Q) shown in Fig. 2b is prominent in l-SiO251, suggesting that Si–O covalent bonds are strong even in the liquid (2019 K). This behavior is consistent with the small differences in density (amorphous: 2.21 g cm-3, liquid: 2.10 g cm-3) and Si–O coordination number, NSi-O (amorphous: 4.0, liquid: 3.9) between them. Considering the melting point of Si (1687 K)52 and SiO2 (1996 K)53, we have noted that the S(Q) of l-Si markedly changes at only 100 K above the melting point, although the temperature-dependent change in the diffraction pattern is small in l-Si32 (see Fig. S1). The marked change in structure above the melting point suggests a large difference in configuration entropy, which is evidence of the poor amorphous-forming ability. On the contrary, the S(Q) of l-SiO2 is similar to that of g-SiO251. We assume that the change in the S(Q) is also an indicator of the amorphous-forming ability.

Figure 2
figure 2

(a) X-ray structure factors, S(Q), for a-Si8 and l-Si (1770 K)29.(b) X-ray structure factors, S(Q), f or a-SiO229 and l-SiO2 (2019 K)51, together with neutron S(Q)30 for a-SiO2. Scattering vector Q is scaled by multiplying by QrA-X (first coordination distance) obtained by a Fourier transform of S(Q). rA-X = rSi-Si = 2.35 Å and 2.45 Å for a-Si and l-Si, and rA-X = rSi-O = 1.61 Å and 1.62 Å for a-SiO2 and l-SiO2, respectively.

Recently, Shi and Tanaka have reported a comparison of S(Q) of several disordered materials with a tetrahedral motif25, by scaling rA-X and rA-A (cation–cation distance). We also display the S(Q) of Si and SiO2 in Fig. 3. As can be seen in the figures, a-Si exhibits an FSDP when the S(Q) is scaled by rA-A, but not when scaled by rA-X. The scaling of Q is still an open question, but, recently Salmon and Zeidler have reported that a-Si does not have an FSDP12.

Figure 3
figure 3

Structure factors S(Q) for a-Si and a-SiO2 scaled by a rA-X (a) and (b) rA-A29. rA-X for a-SiO2 is 1.61 Å, and rA-A for a-SiO2 and a-Si are 3.10 Å and 2.35 Å, respectively.

To understand the relationship between the diffraction peaks and the topology, we calculated ring size distributions on the basis of King and primitive criteria. The primitive and King ring size distributions for crystalline (c)-Si and SiO2 (cristobalite) calculated from crystal structures, and for structure models for liquid and amorphous Si and SiO2 obtained by computer simulations, are shown in Fig. 4 and Fig. S2, respectively. Note that the ring criterion does not affect the overall shape of the distributions. The crystalline phases of both Si and SiO2 exhibit only sixfold rings, in which the ring comprises six silicon atoms in c-Si, whereas it comprises six silicon atoms and six oxygen atoms (AXAX rings) in c-SiO2. The ring size distributions of a-Si and a-SiO2 are broad, and the sixfold rings are dominant for both a-Si and a-SiO2. The ring size distributions of a-Si are consistent with the result of a previous study reported by Opletal et al.54. It is stressed that the ring size distributions of l-Si are significantly different from those of a-Si. The large fraction of threefold rings in the liquid is due to the increased coordination number associated with a significant increase in the density of the liquid. This behavior of the ring size distribution is in line with the difference in diffraction data shown in Fig. 2. On the other hand, the ring size distributions of l-SiO2 (2373 K) are identical to those of a-SiO2, which is consistent with the similarity of the density and coordination number between them. The broad ring size distributions of a-SiO2 and l-SiO2 are a signature of topological disorder according to Gupta and Cooper55, and we previously hypothesized that the width of the ring size distribution is an indicator of the amorphous-forming ability28,56 when short-range structures are the same. Indeed, it appears that a-Si exhibits a narrower ring size distribution than a-SiO2, demonstrating that topological order/disorder is a suitable descriptor for understanding the amorphous-forming ability.

Figure 4
figure 4

taken from ref. 29. It is confirmed that the ring size distribution of l-Si obtained by RMC modeling is identical to the result obtained by density functional (DF)–MD simulation shown in Fig. S3. Since the Si–O-Si bond angle is nearly straight2 in SiO2, we define that a sixfold ring consists of six SiO4 tetrahera, which allows us to decorate oxygen atoms in the analysis. Note that it is not appropriate to compare the ring statistics of l-Si with those of a-Si and SiO2, since the structure of l-Si is beyond the scope of the corner-sharing tetrahedra motif owing to a significantly increased Si–Si coordination number.

Primitive ring size distributions in a series of Si and SiO2. The results of a-Si, l-Si (1770 K), and l-SiO2 (2373 K) are calculated using a structural model obtained by combined molecular dynamics (MD)–reverse Monte Carlo (RMC) modeling29, RMC modeling29, and MD simulation31, respectively. The results of a-SiO2 are

The formation of rings is an important structural feature in covalent amorphous materials; hence, it results in the generation of cavities (empty spaces). We visualize a cavity (highlighted in green) of a-SiO2 together with the cavity size distributions in a-SiO2 and l-SiO2 (2373 K) in Fig. 5. The cavity volume ratio of a-SiO2 is approximately 33%29. The large cavity volume ratio is presumed to be a signature of covalent amorphous materials. Note that the cavity volume ratio of l-SiO2 is comparable to that of a-SiO2 owing to the small difference in density between them. The cavity size distributions in SiO2 indicate that both a- and l-SiO2 have a large cavity in disordered structure. Figure 6 visualizes cavities in a-Si and l-Si together with the cavity size distributions. The cavity volume ratio of 22% for a-Si is much smaller than that of 33% for a-SiO2 shown in Fig. 5, but we observe a significant reduction in cavity volume ratio in l-Si (from 22% in a-Si to 2% in l-Si), which is associated with the increase in density (2.30 to 2.57 g cm−3) and coordination number NSi-Si (3.9 to 5.7). It is worth mentioning that the cavity size distributions indicated that a large cavity in a-Si is squeezed and only small cavities are observed in the liquid. This behavior suggests a poor amorphous-forming ability of Si, which is in line with the fact that a-Si can only be formed in a thin film.

Figure 5
figure 5

(a) Visualization of cavities in a-SiO229 and (b) Cavity size distributions in a- and l-SiO2 (2373 K). The cavity volume ratio of l-SiO2 was calculated on the basis of the atomic configuration obtained by MD simulation. Only the a-SiO2 configuration is shown, because it is difficult to distinguish between a-SiO2 and l-SiO2 by visual inspection.

Figure 6
figure 6

(a) Visualization of cavities and (b) cavity size distributions in a-Si and l-Si (1770 K).

Shi and Tanaka have recently focused on the symmetry of SiSi4 tetrahedra in a-Si and SiSi4 hyper-tetr ahedra in a-SiO225. We apply this analysis to crystalline and amorphous Si and SiO2. Figure 7 shows q values of a series of Si and SiO2. Both c-Si and c-SiO2 have perfect tetrahedra as indicated by the average q value of 1. A comparison between a-Si and a-SiO2 demonstrates that SiSi4 hyper-tetrahedra in a-SiO2 are highly distorted in comparison with regular SiSi4 tetrahedra in a-Si, as revealed by the average q values of a-SiO2 being much smaller than that for a-Si of 0.90. In addition, the distribution of q values is broader in a-SiO2. It is also suggested that hyper-tetrahedra in l-SiO2 are slightly distorted in comparison with those in a-SiO2, whereas the average q value of l-Si is the smallest owing to the formation of fivefold and sixfold Si in the liquid. The profile of the q-parameter for l-Si is consistent with the profile reported in a previous FPMD study49,50. Therefore, the structural motif in the liquid phase is completely different between Si and SiO2 as indicated by the behavior of the density, since the density difference between a-SiO2 and l-SiO2 is small, but large between a-Si and l-Si.

Figure 7
figure 7

SiSi4 tetrahedral order parameter q for a series of silicon and silica. Average q values are given in round brackets. The SiSi4 q parameter of a-SiO2 remains almost the same in densified glass2. The SiO4 q parameters for c-, a-, and l-SiO2 are shown in Fig. S4. It is confirmed that the symmetry of SiO4 in a-SiO2 is better than that in a-Si, and SiO4 in l-SiO2 is more distorted.

Thus, we have revealed the differences in terms of diffraction peaks and the topology by a combination of diffraction measurements and topological analyses. We show that the diffraction data and the topology for a-Si and l-Si are very different. This result is in sharp contrast to SiO2, in which both the diffraction pattern and topology are identical for a-SiO2 and l-SiO2. We are confident that this behavior is the reason why SiO2 is a glass former and Si is not. In other words, the tetrahedral corner-sharing network of AX2, in which A is a fourfold cation and X is a twofold anion, as indicated by the FSDP, is an important motif for the amorphous-forming ability that can rule out a-Si as a good amorphous former. This concept is consistent with the fact that it is impossible to form an elemental bulk amorphous material using melt quenching technique57.

As mentioned in the introduction, the FSDP and PP are suitable descriptors of an amorphous network. The typical behavior of the FSDP and PP under high pressure shown in Fig. 1, in which the FSDP diminishes with a peak shift to a higher Q and the PP becomes sharp with increasing pressure at room temperature, without any peak shifts (Fig. 1a). Onodera et al. have recently reported the evolution of the FSDP with a peak shift in a-SiO2 owing to an LDA-HDA transition induced by heating under high pressure (see Fig. 1b)2. No such significant peak shift can be observed in a-Si under pressures owing to the very strong fully tetrahedral covalent network associated with the absence of an FSDP (note, however, that a first-order-like transformation to high-density a-Si occurs above 10 GPa58,59,60,61). Therefore, it is suggested that the first peak in the S(Q) of a-Si can be assigne d to a PP, because the peak does not show any shifts under high pressures. This interpretation supports that FSDP is a signature of amorphous formin g ability in amorphous network system, although FSDP is not a signature of network29. On the other hand, the FSDP of liquid phosphorus diminishes under high pressures and temperatures, which is associated with the network formation of P4 molecules by transition62,63.

In this article, we demonstrate that topological analyses are powerful tools for understanding the diffraction peaks of disordered materials. Our conclusion fully supports the discussion of Price et al. in 198813. Understanding the diffraction peaks of disordered materials via a series of topological analyses will give rise to the capability to forge a new path for designing novel functional disordered materials.